Tarea Visualización Federico Daverio

¿El salario mínimo impacta negativamente sobre los niveles de empleo? ¿A mayor costo del trabajo, las leyes de la demanda y oferta neoclásica determinan una caída en los niveles requerido de esto y un aumento del desempleo?

Marco Teórico

Según la teoría neoclásica, las fuerzas de oferta y demandan que definen el salario son afectadas cuando se intenta forzar el nivel de retribución. Por otro lado, a pesar de que hay personas que pueden ser beneficiadas con un sueldo mayor, los desplazados son aquellos de productividad baja, es decir, personas recién egresada o con poca experiencia, adultos con baja escolaridad. Por lo tanto, hay consecuencias sociales importantes desde un punto de vista social en la definición del salario mínimo.

Según la evidencia empírica los efectos de un aumento en el salario mínimo serán positivos en la medida que los incrementos sean prudentes y que el salario mínimo inicial sea bajo. En particular, son las condiciones de la economía mexicana.

Objectivo

El objetivo de este estudio es analizar, a través de un método cuasi experimental, el impacto que tuvo en los niveles de empleo el aumento de salario mínimo (SM) implementado en los primeros meses del 2019. En este entonces el SM de 43 municipios del norte del país pasó de 88 pesos al día a 176$. Contemporáneamente, en esta misma área geográfica, se recortó el impuesto sobre el valor añadido que pasó del 16% al 8%. Esta política laboral y retributiva fue circunscrita a la zona norte del país ofreciéndonos las condiciones perfectas para el análisis de un experimento natural.

En particular para evaluar el impacto del conjunto de políticas implementadas construiremos un control sintético para cada municipio tratado considerando como donadores los demás municipios mexicanos. Esto nos permitirá determinar un contrafactual que describa la hipotética trayectoria del empleo en ausencia del cambio en la política de salario mínimo. La heterogeneidad descrita se aprovechará asignando los municipios del norte al grupo tratado y los demás al de control. El delta entre el control sintético y la trayectoria real del empleo en los municipios tratados representará el impacto de la política publica bajo los supuestos metodológicos ilustrados.

Un estudio del 2012 de Campos et al.^[1]aprovechó otra heterogeneidad para evaluar el impacto del salario mínimo en el empleo por medio de Diferencias en Diferencias. Las principales conclusiones a que llegaron fueron que aumentaron los salarios reales, especialmente de los de menores ingresos, el impacto de la política no tuvo efectos significativos en la transición de empleo a desempleo, se redujo la probabilidad de los trabajadores formales hacia la informalidad y finalmente no encontraron efectos sobre los trabajadores formales.

Aprendizaje computacional

¿Como nos puede ayudar el aprendizaje computacional para la realización de este estudio? Hay dos aspectos principales donde resulta fundamental el utilizo de los algoritmos de machine learning para poder efectuar el estudio. El primero es la construcción de cluster de municipios homogéneos adentro de los cuales aplicar el método de control sintético. Esto nos permitirá de disminuir la cardinalidad del conjunto de los donadores permitiéndonos un resultado más preciso y robusto (Abadie et al. 2007)^[2]. Además, siendo que el método de control sintético ocupa una doble optimización anidad que Python resuelve por medio de métodos cuasi-newtonianos de aproximación numérica se disminuyen notablemente las exigencias computacionales.

El segundo ámbito donde el machine learning puede ayudarnos es disminuir la dimensionalidad respecto a las variables explicativas o por lo menos darnos un “insight” sobre las variables observables más importantes para explicar nuestra variable dependiente, en nuestro caso el nivel de empleo. Esto nos permitiría de nuevo ventajas computacionales y una mayor comprensión de las correlaciones entre las observables y la predicha, aunque podría complicar la intuición económica del modelo. Finalmente podríamos utilizar técnicas de aprendizaje computacional también para inicializar la matriz de pesos para las características observables, así de tomar una decisión data-driven respecto a las variables explicativas más relevante para el estudio de nuestro experimento natural.

Visualización

En esta primera fase utilizaremos técnicas de visualización para efectuar un primer acercamiento al análisis de la base de datos. Utilizaremos análisis factorial y de componentes principales para estudiar la relación y correlación entre la variable dependiente y las observables para entender si hay patrones o elementos que resultan más explicativos de otros.

La base de datos

La base de datos empleada fue obtenida de microdatos del IMSS que sucesivamente fueron agrupados a nivel municipio para cuestiones de privacidad y simplificar la intuición. En particular la base de datos contiene los siguientes campos de datos:

  • Y variable dependiente: logaritmo del nivel de empleo a nivel municipal (datos del IMSS)
  • X variables observables/independientes:
    • 9 de porcentaje de trabajo en cada sector
    • 7 de tamaño del empleado
    • 2 porcentaje de trabajo por género
    • 13 porcentaje de trabajo por edad
    • 6 porcentaje de trabajo por nivel de salario
    • 1 poblacional
    • 1 índice de servicio
    • 2 porcentaje de población por edad
    • 1 porcentaje de población con educación básica
    • 1 porcentaje de población con título de graduación

En total tenemos 39 predictores y 1 variable dependiente que utilizaremos en nuestro análisis. En particular para el logaritmo del empleo tenemos su valor para cada municipio para 70 trimestres y en particular para 10 que son sucesivos a la implementación de la política salarial.

Dejaremos todo el código que generó el DB, que fue construido agrupando los datos micro y agrupado y sucesivamente acomodados según las necesidades de los algoritmos de optimización del control sintético, y las visualizaciones al fondo del texto para no sobrecargar la exposición.

Exploración visual

Como variable dependiente estaremos utilizando el logaritmo del nivel de empleo, esto nos permitirá obviar parcialmente a eventuales problemas de heteroscedasticidad. Aprovechando las series de tiempo que disponemos, nos enfocaremos inicialmente en el análisis de la trayectoria del empleo en el tiempo para las unidades tratadas, ósea los municipios donde hubo un incremento en el salario mínimo.

In [2]:
from IPython.display import Image
Image(filename='Figure_1.png')
Out[2]:

En esta primera grafica hemos ploteado el andamiento del logaritmo del empleo en el lapso de tiempo considerado por las unidades tratadas. Podemos notar como no hay tendencias comunes ni en el periodo pretratamiento ni en el postratamiento. Podemos por lo tanto suponer que hay diferencias estructurales marcadas en los municipios analizados que no se pueden suponer heterogéneos. A pesar de esto podemos notar que en ninguno de los municipios interesados por la política publica hubo caídas repentinas en los niveles de empleo después que fue aumentado el salario mínimo.

Para mejorar nuestra comprensión de las diferencias en la evolución del empleo en los municipios en examen graficamos cada trayectoria por separado.

In [3]:
from IPython.display import Image
Image(filename='Figure_2.png')
Out[3]:

Notamos de nuevo que no hay trend evidentemente negativo en la mayoría de los municipios después de la adopción del nuevo salario mínimo. En algunos municipios la caída del empleo pareciera portarse ya desde trimestres anteriores. Solo en el municipio 35 y 23 pareciera haber una disminución del empleo después de t=60 aunque podría depender de otros factores coyunturales porque vemos que esta trayectoria se asemeja a la de los municipios 23 que tuvo un andamiento similar, pero en un punto del tiempo donde no hubo cambios en la política salarial. En la primera grafica vamos a dejar el eje de las y libre para poder apreciar mejor la trayectoria de cada municipio.

In [5]:
from IPython.display import Image
Image(filename='Figure_3.png')
Out[5]:

En el segundo vamos hemos utilizado la misma escala en el eje vertical para poder apreciar la diferencia en magnitud y tener una mejor comparación entre municipios. En apéndice se dejan los gráficos relativos a estas trayectorias, pero sin la aplicación de la función logarítmica a la variable dependiente.

Ahora bien, queremos analizar si hay diferencias no idiosincráticas y estructurales entre los municipios que hacen parte del grupo tratados, donde hubo el aumento del salario mínimo, y los de control. Para hacer esto analizamos las variables observables que hacen parte de nuestra base de datos. Hemos realizado un histograma que refleja la distribución en frecuencia para cada una de estas. En color azul podemos ver las distribuciones para el grupo de control mientras en naranjas las del grupo tratado. Recordamos que tenemos 41 municipios donde se implementó la política salarial y 1753 que se mantuvieron en el estatus quo.

In [9]:
from IPython.display import Image
Image(filename='Figure_6.png')
Out[9]:
In [10]:
from IPython.display import Image
Image(filename='Figure_7.png')
Out[10]:

Podemos notar que efectivamente parecen haber diferencias estructurales entre los municipios tratados y los de control. Esto se debe sobre todo a la grande heterogeneidad geográfica social y económica de México. La frontera norte, a la cuya franja pertenecen todos los municipios tratados, tiene características peculiares que difieren del resto del país. Además, vimos que también comparando entre los mismos municipios del norte había heterogeneidad en la variable dependiente, esto se debe a la granularidad de los datos que lleva casi inevitablemente a grandes varianzas y diferencias antes del agrupamiento. Estas observaciones justificarán el utilizo de un algoritmo de clustering en cuanto nos permite encontrar subconjuntos de municipios (tratados o no) más homogéneos con respecto a las variables explicativas consideradas.

En apéndice podemos encontrar los relativos histogramas de densidad (density=T que confirma los insights obtenidos con la distribución de frecuencia simple.

Ahora bien, graficamos las medias por cada categoría tanto para el grupo tratado como el grupo de control.

In [18]:
from IPython.display import Image
Image(filename='Figure_10.png')
Out[18]:

Podemos apreciar de nuevo que hay diferencias significativas entre los dos grupos y además hay magnitudes/escalas distintas entre las distintas variables observables, en particular por lo que concierne la 33 y la 34. Por este motivo hemos procedido a una estandarización de los datos con la función dedicada de sklearn.

Una vez efectuada la estandarización exploramos la posibilidad de disminuir la dimensión en X siendo que esto nos daría ventajas desde un punto de vista computacional y nos ayudaría a una mejor comprensión de la relación entre olas distintas variables observables y la dependiente. En particular para hacer esto efectuaremos una categorización de la variable dependiente definiendo unos niveles de empleo considerando todos los municipios de la muestra:

  • empleo alto (1): agrupa un tercio de los municipios que tienen los niveles de empleo relativamente más altos dada la muestra (tercer tercíl de la distribución de la variable dependiente)

  • empleo alto (2): agrupa un tercio de los municipios caracterizados por niveles de empleo “medianos” (segundo tercíl de la distribución de la variable dependiente)

  • empleo bajo (3): agrupa un tercio de los municipios que tienen niveles de empleos más bajos con respecto a los demás municipios (primer tercíl de la distribución de la variable dependiente)

Al fin de etiquetar los municipios con respecto a sus niveles de empleo hemos considerado los valores promedios del empleo de las distintas unidades de análisis por todos los periodo pretratamiento (t=60).

El primer análisis realizado fue un análisis factorial. El análisis factorial es un método estadístico que se utiliza para describir la variabilidad entre las variables correlacionadas observadas en términos de un número potencialmente menor de variables no observadas llamadas factores. Lo que buscamos es que las variaciones en las 39 variables observadas se reflejen principalmente en las variaciones de dos o tres variables no observadas (subyacentes/latentes). El análisis factorial busca tales variaciones conjuntas en respuesta a variables latentes no observadas. Las variables observadas se modelan como combinaciones lineales de los factores potenciales, más los términos de "error". Queremos por lo tanto identificar estas variables latentes y encontrar visualmente unas homogeneidades entre los municipios que pertenecen a la misma categoría con respecto a las variables subyacentes identificadas.

Empezamos realizando el análisis factorial con base a dos variables latentes potenciales.

In [19]:
from IPython.display import Image
Image(filename='Figure_11.png')
Out[19]:

Podemos notar como no haya una distribución sistemática distintas de las diferentes categorías en el plan cartesiano. Esto nos induce a pensar que la reducción a solo 2 dimensiones no sea suficiente para explicar la correlación entre las observables de manera tal que sea significativa en la definición de las tres categorías definidas por la variable dependiente. Volvemos a efectuar el análisis factorial considerando tres variables latentes buscando tener mejores “insights”.

In [20]:
from IPython.display import Image
Image(filename='Figure_12.png')
Out[20]:
In [24]:
from IPython.display import Image
Image(filename='Figure_12B.png')
Out[24]:

Aun viendo el espacio tridimensional desde dos perspectivas resulta dificultoso la identificación de una clara distinción en la distribución de los municipios respecto a los distintos niveles de empleo. Esto se podría dar en parte por la cantidad de observaciones y las pocas categorías identificadas y del otro por la varianza entre las observable. A pesar de esto parece haber una mejor división entre los municipios con distintos niveles de empleo, viendo las dos perspectivas notamos que desde un lado se ve una agrupación homogénea de municipios con altos niveles de empleo y del otro resalta el subconjunto de municipios con pocos empleos. Resta, pero complicado entenderlo por medio de la gráfica 3D.

Esto respondería a los supuestos teóricos que prevé que al aumentar del número de variables latentes podemos obtener una mejor aproximación de la varianza original muestral y una mejor explicación de la variable categórica.

Intentamos ahora utilizar otra técnica de reducción dimensional, el análisis de componentes principales (PCA) para ver si podemos obtener resultados mejores. La PCA es una de las técnicas de aprendizaje no supervisado, las cuales suelen aplicarse como parte del análisis exploratorio de los datos permitiéndonos reducir la dimensionalidad, perdiendo la menor cantidad de información (varianza) posible: cuando contamos con un gran número de variables cuantitativas posiblemente correlacionadas (indicativo de existencia de información redundante), la PCA permite reducirlas a un número menor de variables transformadas que expliquen gran parte de la variabilidad en los datos. Hay diferencias sustanciales teóricas de este método respecto al análisis factorial^[3] y la elección entre los dos métodos tendría ser hecha a priori. Por fines didácticos e ilustrativos en nuestro caso aplicaremos ambas.

Empezaremos considerando una reducción a dos componentes.

In [23]:
from IPython.display import Image
Image(filename='Figure_13.png')
Out[23]:

De nuevo podemos notar como la reducción a solo 2 dimensiones no sea suficientes para describir el conjunto de datos en termino de nuevas variables no correlacionada de manera satisfactoria. Por lo tanto, tampoco la descomposición por componente principales parece ayudarnos en la visualización de los datos. De hecho, observamos que no hay subdivisiones claras en el espacio bidimensional entre los municipios de distintas categorías de empleo. Probamos por lo tanto a utilizar 3 componentes principales para explicar la varianza muestral.

In [25]:
from IPython.display import Image
Image(filename='Figure_14.png')
Out[25]:
In [26]:
from IPython.display import Image
Image(filename='Figure_15.png')
Out[26]:

De nuevo podemos notar que la inspección visual a 3 dimensiones parece obtener mejores resultados. Las dos perspectivas utilizadas para el análisis, misma que se ocuparon precedentemente, ahora nos muestran los municipios con niveles medio y alto de empleo, aunque no parece encontrarse una subdivisión sustancialmente mejor respecto a la que se obtuvo con el análisis de factores. Podemos pensar que las dificultades encontradas por los algoritmos en la disminución de dimensiones respecto a las variables explicativas se deben a una categorización demasiado estricta en la variable dependiente (el logaritmo del nivel del empleo). Por lo tanto, volvemos a aplicar los algoritmos precedentemente utilizados, pero con 5 categorías:

  • nivel de empleo alto: 1
  • nivel de empleo medio-alto: 2
  • nivel de empleo medio: 3
  • nivel de empleo medio-bajo: 4
  • nivel de empleo bajo: 5

Se efectuaron los análisis de nuevo por componente principales y por análisis factorial. Se reportan las relativas graficas a continuación.

Análisis factorial:

In [53]:
from IPython.display import Image
Image(filename='Figure_16.png')
Out[53]:
In [54]:
from IPython.display import Image
Image(filename='Figure_17.png')
Out[54]:
In [55]:
from IPython.display import Image
Image(filename='Figure_18.png')
Out[55]:

De nuevo notamos que la representación tridimensional, obtenida con 3 dimensiones, parece obtener resultados mejores a simple inspección visual.

Análisis por componente principales:

In [56]:
from IPython.display import Image
Image(filename='Figure_19.png')
Out[56]:
In [57]:
from IPython.display import Image
Image(filename='Figure_20.png')
Out[57]:
In [58]:
from IPython.display import Image
Image(filename='Figure_21.png')
Out[58]:

Finalmente, también el análisis por componentes principales considerando las 5 categorías no parece obtener resultados mejores en la representación bidimensional, y en la tridimensional confirma los resultados obtenidos precedentemente.

En conclusión, la visualización nos permitió un primer acercamiento a las variables que hacen parte de la base de datos, permitiéndonos una mejor comprensión de la misma, aunque no logramos llegar a resultados concluyentes respecto a la disminución de dimensionalidad.

En esta ultima gráfica, que se reporta a continuación, notamos como varían los valores promedio de los municipios tratados y de control después de la estandardización.

In [59]:
from IPython.display import Image
Image(filename='Figure_22.png')
Out[59]:

Se confirma finalmente una diferencia estructural entre los municipios tratados y los de control, esto representará nuestra justificación teórica la utilizo de los algoritmos de clustering para crear subconjuntos homogéneos adentro de los cuales utilizar el método de control sintético.

Apéndice

In [6]:
from IPython.display import Image
Image(filename='Figure_4.png')
Out[6]:
In [ ]:
Fig. 1-A Trayectorias niveles de empleo municipios tratados
In [7]:
from IPython.display import Image
Image(filename='Figure_5.png')
Out[7]:

Fig. 2-A Trayectorias niveles de empleo municipios tratados (eje y compartido)

In [16]:
from IPython.display import Image
Image(filename='Figure_8.png')
Out[16]:
In [17]:
from IPython.display import Image
Image(filename='Figure_9.png')
Out[17]:

Code

In [46]:
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import FactorAnalysis
from sklearn import decomposition
from sklearn.cluster import KMeans, AgglomerativeClustering
from sklearn.metrics import silhouette_score, calinski_harabasz_score, davies_bouldin_score


#el objetivo de este estudio encontrar un control sintentico para evaluar si hubo impactos significativos
#en el nivel de empleo causados por el aumento del salario minimo. Utilizaremos un metodo cuasi experimental 
#aprovechando de una heterogeneididad en la implementacion de los salarios minimos se consideraron los municipios
#de la frontera norte como grupo tratado y los demas como grupo de control

#siendo particularmente intensivo el calculo del control sintetico cuando hay muchas unidades tratadas si quiere
#clusterizar segun caracteristicas observables homogeneas los municipio mexicanos para efectuar controles sinteticos en muestras más pequeñas

#en esta primera fase analizaremos visualmente la posibilidad de crear estos cluster por medio de diferentes tecnicas 
#de reducción dimensional

#importamos los datos
data1 = pd.read_csv('C:\DAVE2\CIDE\Investigacion\Python mio\mun_data_imss_november_v2.txt',sep="\t", header=None)

#vamos a remplazar los datos indicados como '-999999' para nan siendo que así vienen codificados los 
#missing values en el databse del IMSS
data = data1.replace(-999999, np.nan)

data.shape
#nuestro probelmas esta asi definido
#Y = nivel de empleo a nivel municipal (datos del IMSS)
#X = variables observables/indipendientes: 9 de porcentaje de trabajo en cada sector, 7 de tamaño del empleado
#2 porcentaje de trabajo por genero, 13 porcentaje de trabajo por edad, 6 porcentaje de trabajo por nivel de salario
#1 población. 1 indice de servicio, 2 porcetaje de poblacion por edad, porcentaje de poblacion con educacion basica
#porcentaje de poblacion con titulo de graduacion

#en total tenemos 39 predictores

#definimos la estructura de la base de datos

#la base de datos es asi organizada
#tenemos 74 columnas
# 1 columna es el id de los municipios (unidades tratadas)
U=1794 #numero de municipios totales
Ut=41 #numero de municipios tratados. son los municipios de la zona norte donde se subió el salario minimo
T=70 #maximo numero de periodos a disposición
t0=61 #periodo cuando hubo la intervención
Ta0=T-t0 + 1; #nuemero de periodo post intervencion (el periodo cuando hubo el aumento de salario minimo es considerado como post intervencion)
Tb0=T-Ta0 #periodos pre intervencion

#creamos una matrices con las variables observacionales que nos servirán para predecir la variable dependientes
Xmat=pd.concat([data.iloc[:,10:18], data.iloc[:,19:25], data.iloc[:,26], data.iloc[:,28:40], data.iloc[:,41:46],data.iloc[:,53],data.iloc[:,55:57],data.iloc[:,59:61],data.iloc[:,63],data.iloc[:,65]], axis=1)

#creamos una variable que defina el numero de columnas de la nueva matriz
Xn = len(Xmat.columns)

#creamos una matriz tridimensional donde guardaremos los valores de la variables observables
#para el PERIODO PRE TRATAMIENTO de las VARIABLES EXPLICATIVAS para las UNIDADES TRATADAS
X1=np.zeros((Tb0,Xn,Ut),dtype=float)
i=1

while i<=Ut:
    X1[:,:,i-1]=Xmat[((U-Ut)*T)+(T*(i-1)):((U-Ut)*T)+(T*(i-1)+Tb0)]
    i+=1
    
#creamos una matriz bidimensional con los valores PROMEDIOS de las variables observables
#para el periodo PREINTERVENCION para las UNIDADES TRATADAS
#(se deriva de la matriz precedentemente definida como X1)
Z1=np.zeros((Xn,Ut),dtype=float)

i=1
while i<=Ut:
    Z1[:,i-1]=np.nanmean(X1[:,:,i-1], axis=0).T
    i+=1



#creamos una matriz tridimensional donde guardaremos los valores de la variables observables
#para el PERIODO PRE TRATAMIENTO de las VARIABLES EXPLICATIVAS par a las UNIDADES DE CONTROL    
X0=np.zeros((Tb0,Xn,U-Ut),dtype=float)
i=1

while i<=(U-Ut):
    X0[:,:,i-1]=Xmat[T*(i-1):T*(i-1)+Tb0]
    i+=1

#creamos una matriz bidimensional con los valores PROMEDIOS de las variables observables
#para el periodo PREINTERVENCION para las UNIDADES DE CONTROL
#(se deriva de la matriz precedentemente definida como X0)

Z0=np.zeros((Xn,U-Ut),dtype=float)
i=1

while i<=(U-Ut):
    Z0[:,i-1]=np.nanmean(X0[:,:,i-1], axis=0).T
    i+=1    
    
#Ahora vamos a crear las matrices con las 'y' (osea las verables dependientes)
#Valores variable dependiente (el logaritmo del nivel de empleo) para todos los periodos analizados
#Y para los municipio de CONTROL

#en la columna 71 tenemos el logaritmo del nivel de empleo para los distintos municipios
Ymat=data.iloc[:,71].to_frame()

Ymat2=(Ymat).values

Y0=np.zeros((Tb0,U-Ut),dtype=float)
i=1

while i<=U-Ut:
    Y0[:,i-1]=Ymat2[T*(i-1):T*(i-1)+Tb0,0]
    i+=1



#Niveles de empleo pre intervención para las unidades tratadas
Y1=np.zeros((Tb0,Ut),dtype=float)
i=1

while i<= Ut:
    Y1[:,i-1]=Ymat2[(U-Ut)*T+T*(i-1):(U-Ut)*T+T*(i-1)+Tb0,0]
    i+=1

#Niveles de empleo para la ENTERA serie historica para los municipios TRATADOS
E1 = np.zeros((T,Ut),dtype=float)
i=1
while i<= Ut:
    E1[:,i-1]=Ymat2[((U-Ut)*T)+(T*(i-1)):((U-Ut)*T+(T*(i-1)+T)),0]
    i+=1

#Niveles de empleo para la ENTERA serie historica para los municipios DE CONTROL    
E0 = np.zeros((T,U-Ut),dtype=float)
i=1
while i<= (U-Ut):
    E0[:,i-1]=Ymat2[T*(i-1):T*(i-1)+T,0]
    i+=1

#Creamos un vector de años para poder ver la evolución del empleo en los años
years=list(data[1].unique())

#graficamos la evolución del logaritmo del empleo en los años para los municipios tratados
plt.plot(years,E1)
plt.axvline(x=60,linestyle='--')

#creamos ahora una matriz con los niveles promedio de (log) empleo para las unidades tratadas y la de control
E1M=np.nanmean(E1[0:60][:],axis=0)
E0M=np.nanmean(E0[0:60][:],axis=0)

#junatmos las variables dependientes en un unica matriz
EM=np.concatenate((E1M,E0M),axis=0).T
EM= np.reshape(EM, (-1,1)).T

#creamos intervalos para clasificar el nivel de empleo de las distintas municipalidades
pas = (np.amax(EM)-np.amin(EM))/3
minv=np.amin(EM)

EMC=[]

i=0
while i< EM.shape[1]:
    if EM[0][i]<=(minv+pas):
        EMC.append('bajo')
    elif EM[0][i]<=(minv+ 2*pas):
        EMC.append('medio')
    else:
        EMC.append('alto')
    i+=1

my_dict = {i:EMC.count(i) for i in EMC}
print (my_dict)     

#vemos que no es una buena estrategia siendo que tenemos algun outlier esta sesgando hacia abajo los lavores promedios

#ahora bien dividimos los datos por tericles:

q1=np.quantile(EM,0.33)
q2=np.quantile(EM,0.66)

EMC=[]

i=0
while i< EM.shape[1]:
    if EM[0][i]<=q1:
        EMC.append('bajo')
    elif EM[0][i]<=q2:
        EMC.append('medio')
    else:
        EMC.append('alto')
    i+=1

my_dict = {i:EMC.count(i) for i in EMC}
print (my_dict)     

#Ahora si, por construcción tenemos las 3 categorias equidistribuidas en terminos de numero de observaciones

#para nuestra fianlidad codificaremos las cateogiras con 1=alto, 2 = medio, 3= bajo


EMC=np.where(EMC=='medio',2,EMC)
EMC=np.where(EMC=='alto',1,EMC)
EMC=np.where(EMC=='bajo',3,EMC)
EMC=np.where(EMC=='medio',2,EMC)

#el analisis del empleo es hecha con el logaritmo en cuanto nos permite eliminar parcialmente la heterosedasticidad
#ploteamos cada municipio para averiguar el andamiento del log del empleo en el tiempo
fig1=pd.DataFrame(E1).plot(subplots=T,layout=(9,6), figsize=(18,30))
for c in fig1:
   for ax in c:
      ax.axvline(x=60, color='black',linestyle='--', linewidth=1)

#vamos a plotear cada municipio en un distinto subplot con el eje y con la misma escala 
#para poder hacer comparaciones
fig2=pd.DataFrame(E1).plot(subplots=T,layout=(9,6), sharey=T, figsize=(15,25))
for c in fig2:
   for ax in c:
      ax.axvline(x=60, color='b',linestyle='--') 

#podemos analizar tambien el andamiento del empleo en el tiempo (sin el logaritmo)
fig3=pd.DataFrame(np.exp(E1)).plot(subplots=T,layout=(9,6), figsize=(18,30))
for c in fig2:
   for ax in c:
      ax.axvline(x=60, color='black',linestyle='--', linewidth=1)

fig4=pd.DataFrame(np.exp(E1)).plot(subplots=T,layout=(9,6), sharey=T, figsize=(15,25))
for c in fig2:
   for ax in c:
      ax.axvline(x=60, color='b',linestyle='--')  

#ahora bien calculamos el nivel promedio del log empleo para todos los municipio pre intervencion
#Creamos un database unico que contendrá los valores promedios de las variables observables de interes
#para todos los municipios tratados y no tratados     
Z = np.concatenate((Z0,Z1), axis = 1)

Z0DF=pd.DataFrame(Z0.T)
Z1DF=pd.DataFrame(Z1.T)
ZDF=pd.DataFrame(Z.T)

#ploteamos las distribución de frecuencias de cada una de las caracteristicas observables de los municipios

Z0DF.hist(layout=(9,6), figsize=(15,25))
Z1DF.hist(layout=(9,6), figsize=(15,25), color='salmon')

#y las relativas densidades

Z0DF.hist(layout=(9,6), figsize=(15,25), density=True)
Z1DF.hist(layout=(9,6),  figsize=(15,25), color='salmon', density = True)

#podemos notar que las distribución de frecuencia entre los municipio tratados y los de control difieren sustancialmente
#podemos por lo tanto pensar que hay diferencia sistematicas y no idiosincraticas


#promediamos las caracteristicas observables no estadardiazadas para los municipios
#del grupo de control y los del grupo tratado
meantrat1 = Z.T[0:40,:].mean(axis=0)
meancontr1 = Z.T[41:,:].mean(axis=0)

#ploteamos 
plt.plot(meancontr1)
plt.plot(meantrat1)
       
#notamos la necesidad de estandardizar dado que las escalas difieren notablemente 
#entre las distintas caracteristicas observable utilizada para el estudio


#estandardizamos las variables observadas para homogeneizar las magnitudes
Clus_dataSet = StandardScaler().fit_transform(Z.T)

#averiguamos que efectivamente hemos estadardizado adecuadamente 
#En el eje 0 deberiamos tener media 0 y varianza 1
Clus_dataSet.var(axis=0)
Clus_dataSet.mean(axis=0)

#Efectivamente es el resultado que obtenemos

Clus_dataSet.var(axis=1)
Clus_dataSet.mean(axis=1)
{'alto': 1791, 'medio': 1, 'bajo': 2}
{'medio': 592, 'alto': 610, 'bajo': 592}
Out[46]:
array([0.16081819, 0.18536073, 0.05403275, ..., 0.16027984, 0.13357847,
       0.01282136])
In [44]:
#Entrenamos el algoritmo para la analisis de factores considerando dos dimensiones


k=2
FA=FactorAnalysis(n_components=k).fit(Clus_dataSet)

#ahora predecimos
pre = FA.transform(Clus_dataSet)

for kl in np.unique(EMC):
    m = EMC == kl
    plt.scatter(pre[m, 0], pre[m, 1])
plt.legend(np.unique(EMC),prop={'size': 5},loc=(1.1,0.0))  

#notamos como la reducción a solo 2 dimensiones no sea suficientes para explicar la correlación entre las observables
#y la variable dependiente

#probamos con tres dimensiones:
k=3
FA=FactorAnalysis(n_components=k).fit(Clus_dataSet)
pre = FA.transform(Clus_dataSet)

fig = plt.figure()
ax = plt.axes(projection='3d')

for kl in np.unique(EMC):
    m = EMC == kl
    ax.scatter3D(pre[m, 0], pre[m, 1], pre[m, 2])
ax.view_init(20, 230)
plt.legend(np.unique(EMC),prop={'size': 5},loc=(1.1,0.0)) 


fig = plt.figure()
ax = plt.axes(projection='3d')


for kl in np.unique(EMC):
    m = EMC == kl
    ax.scatter3D(pre[m, 0], pre[m, 1], pre[m, 2])
ax.view_init(20, 120)
plt.legend(np.unique(EMC),prop={'size': 5},loc=(1.1,0.0))   

#parece haber una mejor división aunque es complicado entenderlo por medio de la grafica 3D
#a pesar de esto podemos pensar que al aumentar el numero de componente podriamos obtener una divsion mejor
Out[44]:
<matplotlib.legend.Legend at 0x1f291be4088>
In [37]:
#intentamos el mismo procedimiento pero ahora por componentes principales:

#Entrenamos el algoritmo para el analisis de componente principales considerando dos dimensiones

k=2
PCA=decomposition.PCA(n_components=k).fit(Clus_dataSet)

#ahora predecimos
pre2 = PCA.transform(Clus_dataSet)

for kl in np.unique(EMC):
    m = EMC == kl
    plt.scatter(pre2[m, 0], pre2[m, 1])
plt.legend(np.unique(EMC),prop={'size': 5},loc=(1.1,0.0))  

#de nuevo podemos notar como la reducción a solo 2 dimensiones no sea suficientes para describir el conjunto
#de datos en termino de nyuevas variables no correlacionada. Por lo tanto tampoco la descomposición por componente principales 
#nos ayuda a tener un mejor insight del problema analizado

#probamos con tres dimensiones:

k=3
PCA=decomposition.PCA(n_components=k).fit(Clus_dataSet)
pre2 = PCA.transform(Clus_dataSet)

fig = plt.figure()
ax = plt.axes(projection='3d')

for kl in np.unique(EMC3):
    m = EMC == kl
    ax.scatter3D(pre2[m, 0], pre2[m, 1], pre2[m, 2])
ax.view_init(20, 230)
plt.legend(np.unique(EMC),prop={'size': 5},loc=(1.1,0.0)) 

#y con un distinto punto de observación

fig = plt.figure()
ax = plt.axes(projection='3d')


for kl in np.unique(EMC3):
    m = EMC == kl
    ax.scatter3D(pre2[m, 0], pre2[m, 1], pre2[m, 2])
ax.view_init(20, 120)
plt.legend(np.unique(EMC),prop={'size': 5},loc=(1.1,0.0))   

#de nuevo parece que obtenemos resultados mejores aunque no podemos asegurar esto
#con base a la visualización
Out[37]:
<matplotlib.legend.Legend at 0x1f2912abfc8>
In [39]:
#averiguamos si por medio de una categorización más precisa de la variable dependiente podemos obtener mejores resultados

q1=np.quantile(EM,0.2)
q2=np.quantile(EM,0.4)
q3=np.quantile(EM,0.6)
q4=np.quantile(EM,0.8)

EMC2=[]

i=0
while i< EM.shape[1]:
    if EM[0][i]<=q1:
        EMC2.append(5)
    elif EM[0][i]<=q2:
        EMC2.append(4)
    elif EM[0][i]<=q3:
        EMC2.append(3)
    elif EM[0][i]<=q4:
        EMC2.append(2)
    else:
        EMC2.append(1)
    i+=1

my_dict = {i:EMC2.count(i) for i in EMC2}
print (my_dict)     

#volvemos a efectuar el analisis por visualizacion

#Entrenamos el algoritmo para la analisis de factores considerando dos dimensiones y las 5 categorias

k=2
FA=FactorAnalysis(n_components=k).fit(Clus_dataSet)

#ahora predecimos
pre = FA.transform(Clus_dataSet)

for kl in np.unique(EMC2):
    m = EMC2 == kl
    plt.scatter(pre[m, 0], pre[m, 1])
plt.legend(np.unique(EMC2),prop={'size': 5},loc=(1.1,0.0))

#notamos como la reducción a solo 2 dimensiones de nuevo no sea suficientes para explicar la correlación entre las observables
#y la variable dependiente
{3: 358, 2: 359, 4: 359, 1: 359, 5: 359}
Out[39]:
<matplotlib.legend.Legend at 0x1f2915e4c08>
In [41]:
#probamos con tres dimensiones:
k=3
FA=FactorAnalysis(n_components=k).fit(Clus_dataSet)
pre = FA.transform(Clus_dataSet)

fig = plt.figure()
ax = plt.axes(projection='3d')

for kl in np.unique(EMC2):
    m = EMC2 == kl
    ax.scatter3D(pre[m, 0], pre[m, 1], pre[m, 2])
ax.view_init(20, 230)
plt.legend(np.unique(EMC2),prop={'size': 5},loc=(1.1,0.0)) 


fig = plt.figure()
ax = plt.axes(projection='3d')


for kl in np.unique(EMC2):
    m = EMC2 == kl
    ax.scatter3D(pre[m, 0], pre[m, 1], pre[m, 2])
ax.view_init(20, 120)
plt.legend(np.unique(EMC2),prop={'size': 5},loc=(1.1,0.0))   

#de nuevo parece haber una mejor división aunque es complicado entenderlo por medio de la grafica 3D
#a pesar de esto de nuevo desde el punto de vista teorica deberíamos haber un mejoramiento en la precisión
Out[41]:
<matplotlib.legend.Legend at 0x1f2917e1988>
In [42]:
#de nuevo intentamos el mismo procedimiento pero ahora por componentes principales:

#Entrenamos el algoritmo para el analisis de componente principales considerando dos dimensiones

k=2
PCA=decomposition.PCA(n_components=k).fit(Clus_dataSet)

#ahora predecimos
pre2 = PCA.transform(Clus_dataSet)

for kl in np.unique(EMC2):
    m = EMC2 == kl
    plt.scatter(pre2[m, 0], pre2[m, 1])
plt.legend(np.unique(EMC2),prop={'size': 5},loc=(1.1,0.0))  

#se confirma la dificultad en encontrar patterns a partir de la visualización

#probamos con tres dimensiones:

k=3
PCA=decomposition.PCA(n_components=k).fit(Clus_dataSet)
pre2 = PCA.transform(Clus_dataSet)

fig = plt.figure()
ax = plt.axes(projection='3d')

for kl in np.unique(EMC2):
    m = EMC2 == kl
    ax.scatter3D(pre2[m, 0], pre2[m, 1], pre2[m, 2])
ax.view_init(20, 230)
plt.legend(np.unique(EMC2),prop={'size': 5},loc=(1.1,0.0)) 

#y con un distinto punto de observación

fig = plt.figure()
ax = plt.axes(projection='3d')


for kl in np.unique(EMC2):
    m = EMC2 == kl
    ax.scatter3D(pre2[m, 0], pre2[m, 1], pre2[m, 2])
ax.view_init(20, 120)
plt.legend(np.unique(EMC2),prop={'size': 5},loc=(1.1,0.0))   
#se confirman los resultados obtenidos con 3 categorias
Out[42]:
<matplotlib.legend.Legend at 0x1f291989708>
In [43]:
#finalmente calculamos el promedio de las caracteristica observables estamdardizadas para los municipios TRATADOS
meantrar=Clus_dataSet[0:40,:].mean(axis=0)

#calculamos el promedio de las caracteristica observables estamdardizadas para los municipios de CONTROL
meantrat=Clus_dataSet[0:40,:].mean(axis=0)

meancontr=Clus_dataSet[41:,:].mean(axis=0)

plt.plot(meancontr)
plt.plot(meantrat)

#Recordamos que hemos etandardizado respecto a todo el grupo analizado no solo el de control por esto 
#cuando consideramos el grupo de control las caracteristicas no tienen media 0
#vemos pero que es congrunte siendo que la varianza es siempre inferior a 1

#finalmente notamos que las caracteristicas del grupo tratado y el grupo de control difieren notablemente
#por casi todas las caracteristicas observables 
Out[43]:
[<matplotlib.lines.Line2D at 0x1f291ac6c88>]

II Parte

In [ ]:
#clustering with k-means, with n_cluster it's possible to modify the number of desired cluster while with n_init you can
#establish the number of loop before the interrumption of the algorithm
#K1=[]
#K=[]
#Kc=[]
#K2 = []
#for k in range (2,30):
    #k_means = KMeans(init = "k-means++", n_clusters = k, n_init = 200).fit(Clus_dataSet)
    #labels = k_means.labels_
    #cl = k_means.predict(Clus_dataSet)
    #cl2 = AgglomerativeClustering(n_clusters=k).fit_predict(Clus_dataSet)
    #_ = calinski_harabasz_score(Clus_dataSet, cl)
    #K1.append(_)
    #n = silhouette_score(Clus_dataSet, cl)
    #K.append(n)
    #Kc.append(calinski_harabasz_score(Clus_dataSet, cl2))
    #n2 = davies_bouldin_score(Clus_dataSet, cl)
    #K2.append(n2)
    


#plt.plot(range(2,30),K1,label='Calinsky k mean')
#plt.plot(range(2,30),K,label='Silhouette k mean')
#plt.plot(range(2,30),Kc,label='Silhouette Aglomerative clustering')
#plt.plot(range(2,30),K2,label='Davies k mean')
#plt.legend()
#assign the label to the different municipality
#Clusdf2 = pd.DataFrame(Clus_dataSet)
#Clusdf2["Clus_km"] = labels
#print(Clusdf2.head(5))




#mean value for all the explicable variable in a group
#print(Clusdf2.groupby('Clus_km').mean())



#add the label to identify the unit treated
#labtrat = pd.DataFrame(np.concatenate((np.array([0]*(U-Ut)),np.array([1]*Ut)), axis = 0)) 
#Clusdf3 = Clusdf2.copy()
#Clusdf3['treat']= labtrat

#resume of the municipality in each group and the unit treated in the same group
#resumen = pd.concat([Clusdf3.groupby('Clus_km').count()['treat'],Clusdf3.groupby('Clus_km').sum()['treat']], axis=1)
#resumen.columns=['# mun in group ', 'mun tr in group']
#print(resumen)

#seleccionamos el primero cluster y hacemos el analisis
#Clusdf3.loc[Clusdf3['Clus_km'] == 0]

Bibliografía

[1]: Raymundo M. Campos, Gerardo Esquivel y Alma S. Santillán (2015). El impacto del salario mínimo en los ingresos y el empleo en México. Cepal- Serie de estudios y perspectivas-México- Nº 162.

[2]: Abadie, Alberto & Diamond, Alexis & Hainmueller, Jens. (2007). Synthetic Control Methods for Comparative Case Studies: Estimating the Effect of California's Tobacco Control Program. Journal of the American Statistical Association. 105. 493-505. 10.1198/jasa.2009.ap08746.

[3]: Santos, Roberta de Oliveira, Gorgulho, Bartira Mendes, Castro, Michelle Alessandra de, Fisberg, Regina Mara, Marchioni, Dirce Maria, & Baltar, Valéria Troncoso. (2019). Principal Component Analysis and Factor Analysis: differences and similarities in Nutritional Epidemiology application. Revista Brasileira de Epidemiologia, 22, e190041. Epub July 29, 2019.https://doi.org/10.1590/1980-549720190041

In [ ]: